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ABSTRACT 

The evolution of halos consisting of weakly self-interacting dark matter 
particles is investigated using a new numerical Monte-Carlo N-body method. 
The halos initially contain kinematically cold, dense r _1 -power-law cores. For 
interaction cross sections a* — a/m p > 10 - 100 cm 2 g" 1 weak self-interaction 
leads to the formation of isothermal, constant density cores within a Hubble 
time as a result of heat transfer into the cold inner regions. This core structure 
is in good agreement with the observations of dark matter rotation curves in 
dwarf galaxies. The isothermal core radii and core densities are a function of 
the halo scale radii and scale masses which depend on the cosmological model. 
Adopting the currently popular ACDM model, the predicted core radii and core 
densities are in good agreement with the observations. For large interaction 
cross sections, massive dark halos with scale radii r s > 1.4 x 10 4 (cm 2 g _1 /a*) 
kpc could experience core collapse during their lifetime, leading to cores with 
singular isothermal density profiles. 

Subject headings: dark matter - galaxies: halos - galaxies: formation - galaxies: 
kinematics and dynamics - methods: numerical 
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1. Introduction 

Cosmological models with a dominant cold dark matter component predict dark 
matter halos with strongly bound, kinematically cold cores (Dubinski & Carlberg 1991, 
Warren et al. 1992, Navarro et al. 1997). Within the core region, the dark matter density 
increases as a power-law p ~ r~ 7 with 7 in the range of 1 to 2 and the velocity dispersion 
a decreases towards the center (Carlberg 1994, Fukushige & Makino 1997). Numerous 
numerical simulations (e.g. Moore et al. 1998, Huss et al. 1999, Jing & Suto 2000), as 
well as analytical theory (Syer & White 1998, Kull 1999), have shown that such a core 
structure follows naturally from collisionless hierarchical merging of cold dark matter halos, 
independent of the adopted cosmological parameters. 

It has recently become clear that on galactic scales the predictions of cold dark matter 
models are not in agreement with several observations. High-resolution calculations by 
Klypin et al. (1999) and Moore et al. (1999) have shown that the predicted number 
and mass distribution of galaxies in galactic clusters is consistent with the observations. 
However, on scales of the Local Group, roughly one thousand dark matter halos should exist 
as separate, self-gravitating objects, whereas less than one hundred galaxies are observed. 
This disagreement can be attributed to the high core densities of satellite dark halos in 
cosmological models which stabilize them against tidal disruption on galactic scales. Mo et 
al. (1998) and lateron Navarro and Steinmetz (2000) found that cold-dark matter models 
reproduce well the I-band Tully-Fisher slope and scatter. They however fail to match 
the zero-point of the Tully-Fisher relation as well as the relation between disk rotation 
speed and angular momentum. Again, this problem can be traced to the excessive central 
concentrations of cold dark halos. Finally, recent observations of dark matter dominated 
rotation curves of dwarf galaxies have indicated shallow dark matter cores which can be 
described by isothermal spheres with finite central densities (Moore 1994, Flores & Primack 
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1994, Burkert 1995, de Blok & Mc Gaugh 1997, Burkert k Silk 1999, Dalcanton & Bernstein 
2000, see however van den Bosch et al. 1999), in contrast to the power-law cusps, expected 
from cosmological models. The disagreement between observations and theory indicates 
that a substantial revision to the cold dark matter scenario might be required which could 
provide valuable insight into the origin and nature of dark matter. 

Motivated by these problems, Spergel & Steinhardt (1999) proposed a model where 
dark matter particles experience weak self-interaction on scales of kpc to Mpc for typical 
galactic densities. They noted that self-interaction could lead to satellite evaporation due 
to the dark particles within the satellites being kicked out by high-velocity encounters with 
dark particles from the surrounding dark halo of the parent galaxy. In order for weak 
interaction to be important on galactic scales, they estimate that the ratio of the collision 
cross section and the particle mass should be of order a wsi /m p fa 1 cm 2 g _1 . 

The Spergel and Steinhardt model has already motivated several follow-up studies. 
For example, Ostriker (1999) demonstrated that weak self-interaction would have the 
interesting side product of naturally growing black holes with masses in the range 10 6 — 10 9 
M in galactic centers. Hogan & Dalcanton (2000) investigated analytically the effect of 
particle self-interactions on the structure and stability of galaxy halos. Moore et al. (2000), 
adopting a gas-dynamical approach, showed that in the limit of infinitely large interaction 
cross sections dark halos would develop singular isothermal density profiles which are not in 
agreement with observations. Mo & Mao (2000) and Firmani et al. (2000) investigated the 
affect of self-interaction on rotation curves. In addition, models of repulsive dark matter 
(Goodman 2000), fluid dark matter (Peebles 2000) and self-interacting warm dark matter 
(Hannestad & Scherrer 2000) have recently been discussed. 

In this paper we will investigate the effect of weak self-interaction on the internal 
structure of cold dark matter halos. If the interaction cross section is not exceptionally 
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large, the dark matter system cannot be treated as a collision dominated, hydrodynamical 
fluid. Section 2 therefore introduces a new numerical Monte-Carlo-N-body (MCN) method 
for weakly interacting particle systems. Initial conditions are discussed in section 3. Using 
the MCN-method, the evolution of weakly self-interacting dark matter halos is investigated 
in section 4. Conclusions follow in section 5. 

2. The Monte-Carlo N-body method 

Within the framework of weak self-interacting, the mean free path A of a dark matter 
particle is determined by A = (per*) -1 , where p is the local dark matter mass density and 
o"* = o- wsi /m p is the ratio between the self-interaction collision cross section a wsi and the 
particle mass m p . If, within a timestep At, the path length / = vAt of a particle with 
velocity v is short compared to A, the probability P for it to interact with another particle 
can be approximated by 

p = l/\ = a *pvAt. (1) 

We use a Monte-Carlo approach in order to include weak self-interaction in a 
collisionless N-body code that utilizes the special purpose hardware GRAPE (GRAvity 
PipE; Sugimoto et al. 1990) in order to determine the gravitational forces between the 
dark matter particles by direct summation. For each particle, a list of its 50 nearest 
neighbors is returned by the boards which allows the determination of the local dark matter 
mass density p. The particle experiences an interaction with its nearest neighbor with a 
probability given by equation (1). Each weak interaction changes the velocities of the two 
interacting particles. Here, due to the lack of a more sophisticated theory, we assume that 
the interaction cross section is isotropic and that the interaction is completely elastic. In 
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this case, the directions of the velocity vectors after the interaction are randomly chosen 
and their absolute values are completely determined by the requirement of energy and 
momentum conservation. 

The computational timestep At must be chosen small enough in order to guarantee 
that the evolution is independent of the numerical parameters, that is the adopted timestep 
and the number of particles. Otherwise, particles with large velocities could penetrate too 
deeply into a dense region like the core of a dark matter halo, violating the requirement 
I << A. Test calculations have shown that At < ^(o^p-u) -1 with rj « 0.1 leads to reliable 
results that are independent of the numerical parameters. 

3. Initial conditions 

Cold dark matter halos form on dynamical timescales. If a* is small enough, the halos 
will achieve an equilibrium state within a few dynamical timescales that is determined by 
collisionless dynamics alone, before self-interaction becomes important. The structure of 
the halos subsequently changes due to self-interaction on longer timescales. This secular 
evolution is similar to the long-term evolution of globular clusters which experience core 
collapse due to gravitational 2-body encounters after virialization. 

We start with an equilibrium model of a virialized dark matter halo and study its 
secular dynamical evolution due to weak self-interaction using the MCN-method. As initial 
condition, a Hernquist halo model (Hernquist 1990) is adopted. Its density distribution 
is p(r) = psj (r/r s (l + r/r s ) 3 ) where p s and r s are the scale density and scale radius, 
respectively. The mass profile is M(r) = Mr 2 /(r s + r) 2 with M the finite total halo 
mass. Assuming hydrostatic equilibrium and an isotropic velocity distribution, the velocity 
dispersion is zero at the center and increases outwards, reaching a maximum at the 
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inversion radius = 0.33r s outside of which it decreases again. A similar structure is seen 
in cosmological simulations (Carlberg 1994, Fukushige & Makino 1997). In general, within 
the interesting region r <r s the Hernquist model provides an excellent fit to the structure 
of cold dark matter halos that result from high-resolution cosmological models. Only in the 
outermost regions do the dark matter halo profiles deviate significantly from the Hernquist 
model, predicting a density distribution that decreases as r~ 3 and a dark halo mass that 
diverges logarithmically (Navarro et al. 1997). Note, that our model neglects any clumpy 
substructure that might exist within dark matter halos (Moore et al. 1999). This should be 
a reasonable approximation for the inner regions where satellites are efficiently disrupted by 
tidal forces. The evolution of weakly self-interacting, clumpy dark halos will be presented 
in a subsequent paper (see also Moore et al. 2000). 

In the following, we will adopt dimensionless units: G=l, r s = 1 and M = 1. The total 
mass and the mean mass density within the inversion radius r\ is Mj = 0.06 and pi = 0.4, 
respectively, leading to a dynamical timescale within of r^ n = 0.8. Most numerical 
calculations have been performed adopting 80000 particles and a gravitational softening 
length of e = 0.002 xr s . Test calculations with 120000 particles did not change the results. 
N-body calculations without weak interaction have shown that the dark halo is stable 
and its density distribution does not change outside of r > 0.006 xr s within 20 dynamical 
timescales. 

4. The evolution of weakly self-interacting dark halos 

Figure 1 shows the evolution of the dark matter density distribution and the velocity 
dispersion profile inside the core region, adopting a collision cross section <r*=10xr^/M s . 
The density distribution initially has the characteristic power-law cusp and the velocity 
dispersion decreases towards the center for r < r\. Within this region, the kinetic 
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temperature inversion leads to heat conduction inwards. The central velocity dispersion 
increases with time and the core expands, resulting in a shallower density distribution. 
After 3 dynamical timescales, an isothermal, constant density core has formed with a radius 
that is of order the initial inversion radius r\. Subsequently, weak interactions between the 
kinematically hotter core and the cooler envelope lead to a flow of kinetic energy outwards 
which causes the isothermal core to contract and heat up further due to its negative specific 
heat, starting a core collapse phase. The calculations are stopped after 16 dynamical 
timescales when the central density and the central velocity dispersion has increased further 
by a factor of 4 and 1.4, respectively. Note that during the core collapse phase, the system 
maintains an isothermal, constant density core with the core radius decreasing with time. 
Overall, the evolution of the dark halo is very similar to the secular evolution of particle 
systems with Hernquist profiles that are affected by gravitational 2-body interactions 
(Heggie et al. 1994, Quinlan 1999). 

Several calculations with different interaction cross sections a* have been performed. 
In all cases, the evolution is similar to that shown in Fig. 1, independent of the adopted 
collision cross section. The timescale r iso for the formation of the isothermal constant 
density core does however depend on er* with 



Tiso w ~^w s - {2) 

In agreement with the calculations of Quinlan (1999) the core collapse timescales are 
roughly an order of magnitude larger than r iso . 

Observations of dark matter dominated dwarf galaxies show a characteristic dark 
matter core structure that can be fitted well by the empirical density distribution (Burkert 
1995) p = po * (r + r ) _1 (r 2 + r^) -1 where po and r are the isothermal core density and 
radius, respectively. Figure 2 compares this profile (solid line) with the core structure of 
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weakly interacting dark halos at t — (dashed line) and after core expansion at t — r iso 
(points with error bars). It is well known that power-law cores do not provide a good fit to 
the observations. An excellent agreement can however be achieved after core expansion if 
one adopts the following core parameters 



r w 0.6r, (3) 
po « 1.54M;r- 3 



where M,- is the initial dark matter mass inside the inversion radius r,- 



5. Conclusions 

The previous MCN calculations have shown that isothermal cores with shallow density 
profiles form naturally in weakly interacting dark halos. The cores have density distributions 
that are in excellent agreement with the observations of dark matter rotation curves in 
dwarf galaxies. The core size is determined by the radius r^ inside which heat is conducted 
inwards, that is where the initial velocity dispersion decreases towards the center. Note, 
that this conclusion should be valid, independent of whether the density diverges as p ~ r _1 
or even steeper (Moore et al. 1998, Jing & Suto 2000) for r <C r s . 

A quantitatively comparison with the observations requires the determination of the 
typical scale parameters r s and M s for dark matter halos. Recent cosmological ACDM 
models (Navarro & Steinmetz 2000) predict that cold dark matter halos with total masses 
M 20 o ~ 10 10 — 10 12 M should have concentrations c = r 200 /r s ~ 20, where M 2 oo is the 
total dark matter mass within the virial radius r 2 oo which denotes the radius inside which 
the averaged overdensity of dark matter is 200 times the critical density of the universe. 
Adopting a Hubble constant h=0.7 leads to r 20 o = 0.02 (Maoo/M©) 1 / 3 kpc « 40 - 200 kpc 
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and with c=20 to scale radii r s w 2 - 10 kpc. For a NFW-profile (Navarro et al. 1997) 
the dark matter mass inside r s is M s 0.1 M 2 oo and the core density is M s /r 3 ks 0.01 
M Q pc~ 3 . In contrast to the Hernquist model with r\ = 0.33r s , the inversion radius of the 
NFW-profiles coincides with the scale radius r\ = r s due to the shallower outer density- 
distribution. According to equation 3, weak interaction in NFW-halos should therefore 
lead to isothermal cores with radii r ~ 0.6 r s w 1.2 - 6 kpc and densities po ~ 1-55 M s 
r~ 3 m 1.5 x 10~ 2 M pc~ 3 . The observations indicate core radii r ~ 2 - 10 kpc with core 
densities p ~ 0.01 M Q pc~ 3 (Burkert 1995), in very good agreement with the theoretical 
predictions. 

In order for dark matter cores to be affected by weak self-interaction, the core expansion 
timescale must be smaller than the age r of the halo: T iso < r ps lOOr^. With equation 
(2) and adopting M s /r 3 ps 0.01M o pc~ 3 , this requirement leads to a minimum value of the 
collision cross section for weak self-interaction to be important 



Note, that this lower limit would be a factor of 25 larger if cosmological models 
underestimate the scale radii of dark halos by a factor of 5. 

Dark matter halos with r iso ps Td yn are likely to have gone through core collapse if 
their ages are r >> r dyn . This condition requires the halo core radii to be larger than 
r s > 1.4 x 10 4 (cm 2 g -1 )/cr* kpc, indicating that more massive halos could have experienced 
core collapse while lower mass halos could still be in the process of core expansion. 

I would like to thank Matthew Bate for providing a subroutine to find nearest 
neighbors using GRAPE and Paul Steinhardt, Ben Moore and Jerry Ostriker for interesting 
discussions and the referee for important comments. Special thanks to David Spergel for 
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pointing out the different dependence of the inversion radius on the scale radius for NFW- 
and Hernquist profiles. 
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Fig. 1. — Figure 1: The evolution of the dark matter density distribution p(r) and of the 
3-dimensional velocity dispersion profile cr(r) in dimensionless units is shown, adopting a 
weak interaction cross section of a* = 10 r 2 s /M s . The upper panels show the phase of core 
expansion with the dotted curves representing the initial conditions and the dashed and 
solid curves corresponding to the evolutionary state of t=l r dyn and t=2.5 r dyn , respectively 
The lower panels show the epoch after core expansion (lower dotted lines: t=2.5 r dyn , upper 
dotted lines: t=5 Td yn ) an d the subsequent phase of core collapse (dashed line: t=10r dyn , 
solid line: t=16 r dyn ). 
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Fig. 2. — Figure 2: The dark matter density distribution as inferred from rotation curves 
of dwarf galaxies (solid line) is compared with halo density distribution at t=0 (dashed line 
) and after core expansion at t = r iso . 



